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Abstract. We develop a path- integral theory to study the angle-resolved photoemission 
spectra (ARPES) of high-T c superconductors based on a two-dimensional model for the 
CuC>2 conduction plane, including both electron-electron (e-e) and electron-phonon (e-ph) 
. interactions. Comparing our result with the experimental one of Bi2Sr2CaCu2 08, we find that 

the experimentally observed isotopic band shift in ARPES is due to the off-diagonal quadratic 
e-ph coupling, whereas the presence of e-e repulsion partially suppresses this effect. 



O . 

Q , 1. Introduction 

The study of high-T c superconductivity is one of the most attractive realms in the last two 
decades. Since the angle-resolved photoemission spectroscopy (ARPES) directly probes the 
electronic occupied states, it has become an important technique to investigate the electronic 
properties of cuprates [T|. Recently, the oxygen isotope effect has been detected with ARPES in 
Bi2Sr2CaCu208 (Bi2212) by two groups[21 [3], and a common feature is noticed that the spectra 
are shifted with the 16 0/ 18 substitution, providing direct evidence for electron-phonon (e-ph) 
coupling in this material. However, since the first report by Gweon et a/. [2], this isotopic band 
shift has become an controversial issue [4], as the observed shift is up to 40 meV, much larger 
than the isotopic energy change of oxygen phonon, ~ 5 meV[5]. Very recently, Douglas et al.[S\ 
repeat the experiment, but find the shift is only 2±3 meV. Thus it turns out to be an interesting 
problem whether the large shift observed by Gweon et al. is possible or not in the cuprates. In 
this paper we shall look into this isotope induced band shift from a theoretical point of view. 

In the Cu02 plane of cuprates, as shown in Fig. 1, the electronic transfer is affected by the 
vibration of oxygen atoms between the initial and final Cu sites, resulting in an off-diagonal 
type e-ph coupling. In order to have an insight into the isotope effect of B42212, we start from a 
half-filled Hamiltonian including the electron-electron (e-e) repulsion and the above mentioned 
off-diagonal e-ph coupling (h = 1 and ks = 1 throughout this paper): 

H = - t(M')(«U, + aJ v %) + ^Evi| + irE + > (!) 

<j,jv 1 <v> V qw ' 

where a] (ai a ) is the creation (annihilation) operator of an electron with spin a at the Cu site I 
on a square lattice (see in Fig. 1). The electrons hop between two nearest neighboring Cu sites, 
denoted by (I, I'), with a transfer energy t(l, I'). U is the strength of Coulomb repulsion between 



two electrons on the same Cu site with opposite spins. The oxygen phonon is assumed to be of 
the Einstein type with a frequency loq and a mass m. A (= 1 + Am/m) is the mass change factor 
of phonon due to the isotope substitution. In the third term, qui is the dimensionless coordinate 
operator of the oxygen phonon locating between the nearest-neighboring Cu sites I and and 
the sum denoted by (I, I') just means a summation over all the phonon sites in the lattice. 
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Figure 1. A Schematic plot of Cu02 conduction 
plane in cuprates. The copper atom (white circle) 
is on the Z-th site of a simple square lattice. The 
oxygen atom (black circle) is located between the 
nearest-neighboring Cu sites, and qu> denotes its 
displacement from the equilibrium position. 



In the conduction plane of Cu02, the electronic hopping integral t(l,V) can be expanded to 
the second order terms with respect to the phonon displacements qu> as 

t(l,l')=t + sq? l ,, (2) 

where to is the bare hopping energy and s the off-diagonal quadratic e-ph coupling constant. 
Here we note the linear e-ph coupling does not occur owing to the lattice symmetry of present 
model. Whereas the inter-site e-e interaction is included in the screened values of to an d s. 



2. Path-integral Monte Carlo 

In this section, we develop a path-integral theory for a model with both e-e and e-ph interactions. 
By making use of the Trotter's decoupling formula, the Boltzmann operator is written as, 



e -f>H 



lim (e~ AH ---e- AH ). (3) 



Applying the Hubbard-Stratonovitch transformation [H] and the Gaussian integral formula[7j, 
we decouple the two-body parts, so that the e-e and e-ph correlated terms are replaced by a 
two-fold summation over the auxiliary spin and lattice configurations, which is the so-called 
path-integral. In this way, the Boltzmann operator is rewritten into the path-integral form as, 

€~ m |pxfr + exp|-^dr[/2(r,x m ,x g )+O(x g )]|n[l^(^^)>(^(^0)i]j , (4) 



h \T, X m , Xq 



[t + sx 2 q (l,l',T) a\ a (T)a Va (T) + a\ l(T (T)ai a {T) 
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^x m (Z,T)[n /T (r) -n u (r)], 
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'dx q (l,l',T)~ 




[ dr \ 



+ ^ xI(1,1',t) 



(5) 
(6) 



Here, x m and x q correspond to the auxiliary spin and lattice field, respectively, and / Vx 
symbolically denotes the integrals over the path x synthesized by x m and x q . A is the time 



interval of the Trotter's formula, = 1/T, and T is the absolute temperature. T + in Eq. (4) is 
the time ordering operator. 

Then we define the time evolution operator [= R(t, x)] as 



R(t, x) = T + exp 



/ dr'h{r [ ',x m ,: 
Jo 



(7) 



In terms of the Boltzmann operator (4) and time evolution operator (7), we define the free 
energy [= <&(#)] of the given path as 



e -/J*(*) = e -Io dTn ^Tr[R((3,x)] . 
While, the partition function (= Z) and total free energy (= $) are given as 

Z = = J Vxe~^ x \ 



(8) 



(9) 



According to Refs. [6] and [7], we also define the one-body Green's function [= G a (lT, I't', x)] 
on a path x as 

G a (lT,lV,x) = -sign(r - T'){T + a la {r)a\ l(T {T')) x , (10) 
where a; CT (r) is the Heisenberg representation of ai a . It is really time-dependent and defined by 

a/ CT (r) = R~ 1 (T,x)ai a R(T,x). (11) 
Meanwhile, the ordinary Green's function [= G a (l,T)] can be obtained by the path-integral as 

G(j (I -l',r-r') = ^j Vxe-^G a (h, I't', x). (12) 

This path-integral is evaluated by the quantum Monte Carlo (QMC) simulation method. 

If the QMC data of Green's function G a (l,T) is obtained, we can immediately calculate its 
Fourier component [= G a (k, r)] as 



1 



(13) 



where k is the momentum of the outgoing photo-electron. From this Fourier component G a (k, r) , 
we derive the spectral function [= Aj(k, oj)] by solving the integral equation 



duj- 3- 

oo 1 + e-P" 



A a (k,oj) 



Finally, the normalized spectral intensity is obtained as, 



I(k,co) 



where the Fermi-Dirac function f(u) = l/[exp(/?u;) + 1] is imposed. 
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Figure 2. (a) The variation of Aer/Aw with s on a 4x4 square lattice, when [7=2.0, 
/3=10, Ao=l, A=2. The filled circles are from QMC, and the solid curve from HFA as a 
guide for eyes, (b) The variations of Aer/Ao; with U on a 4x4 square lattice at (3=10, 
Ao=l and A=2. Three different values of s are used to show the s-dependence of Aer/Ao;. 
The discrete symbols are the results of QMC, and continuous curves by HFA as a reference. 



3. Results and Discussions 

We now present the QMC results on a 4x4 square lattice, where to is set as the unit of 
energy, and ujq=1.0 is used. For the QMC simulation, we impose a little large isotopic mass 
enhancement, Ao=l and A=2, to suppress the numerical error. In this calculation, we determine 
the binding energy ek by the moment analysis of the spectral intensity as ek = / <io;I(k, 
Correspondingly, the isotope induced band shift is calculated by Aek — ek(A) — ek(Ao). 

In Fig. 2(a), we plot the ratio Aer/Ao; versus s, at [7=2.0 and /3=10, where Aer is the 
band shift at the T point of Brillouin zone [kr=(0,0)], and Au> is the isotopic change of phonon 
energy. The filled circles are calculated by QMC, and the solid curve by the mean-field theory 
with Hartree-Fork approximation (HFA) as a guide for eyes. Here both theories figure out an 
increase of Aer/Ao; with s, which means if the e-ph coupling is strong enough, a large band shift 
can be generated in the cost of a small Aw. In Fig. 2(b), the ratio Aer/Ao; versus U are shown 
for three different s's, where the discrete symbols and continuous curves are the QMC and HFA 
results, respectively. It can be seen that the ratio Aer/Ao; increases with s. Meanwhile, for a 
fixed s, the ratio declines slightly as U increases. This behavior indicates that the band shift 
is owing to the e-ph coupling, whereas the presence of U partially reduces this effect. In terms 
of Figs. 2(a) and 2(b), one can see the band shift is actually a measure of the e-ph coupling 
strength in the system. If the result of Ref. [2] is correct, the e-ph coupling must be strong in 
Bi2212. On the contrary, Ref. [3] shows that the coupling cannot be very large. 

4. Conclusion 

In summary, by using the path-integral QMC method, we study the isotopic shift in the ARPES 
of Bi2212 based on a model including both e-e and off-diagonal quadratic e-ph interactions. Our 
calculation demonstrates that the band shift is primarily triggered by the e-ph coupling, while 
the presence of e-e repulsion tends to suppress this effect. 
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